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Prom non-degenerate conducting polymers to dense matter 
in the massive Gross-Neveu model 
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Using results from the theory of non-degenerate conducting polymers like cis-polyacetylene, we 
generalize our previous work on dense baryonic matter and the soliton crystal in the massless Gross- 
Neveu model to finite bare fermion mass. In the large N limit, the exact crystal ground state can 
be constructed analytically, in close analogy to the bipolaron lattice in polymers. These findings are 
contrasted to the standard scenario with homogeneous phases only and a first order phase transition 
at a critical chemical potential. 

PACS numbers: 11.10.Kk,11.10.Lm,11.10.St 



I. INTRODUCTION 

In this paper we continue our investigation of the 
Gross-Neveu (GN) model lj, a 1+1 dimensional, asymp- 
totically free relativistic quantum field theory with La- 
grangian 

C = ^(jrfdp - m )V (i) + \g 2 {^H {i) f. (1) 

The index i runs over N fermion species, and we are only 
concerned with the 't Hooft limit N — ► oo, Ng 2 = const. 
If the bare fermion mass mo vanishes, this model has a 
discrete chiral symmetry ip — > ^ip. During the last few 
years, it has become clear that the phase diagram at fi- 
nite temperature and chemical potential of the massless 
GN model is richer than previously thought Q , the most 
striking novel feature being the appearance of a solitonic 
crystal phase. This phase diagram can be computed ex- 
actly, with many results in closed analytical form 0, Q . 
The corresponding analysis for the massive case (mo 7^ 0) 
has yet to be done. As a first step in this direction, we 
have recently determined the properties of single baryons 
in the massive GN model 0- Here we address the prob- 
lem of dense baryonic matter at zero temperature in the 
massive GN model. 

Interestingly, the GN model has lived a kind of dou- 
ble life in particle physics and condensed matter physics 
over marry years, with only sporadic cross-references (see 
e -g- 0- 13)- On the one hand, it is a soluble model 
for strong interaction physics, exhibiting phenomena like 
asymptotic freedom, chiral symmetry breaking, dynami- 
cal mass generation, meson and baryon bound states . 
On the other hand, it can model quasi-one dimensional 
systems in the vicinity of a half-filled band, in particu- 
lar conducting polymers |Toj . A prominent example is 
iraras-polyacetylene (PA). This polymer (CH) X possesses 
two dimerized, degenerate ground states (Fig. 1) and, 
owing to a number of simplifying assumptions, leads to a 
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FIG. I: Two degenerate, dimerized ground states of trans-PA, 
described theoretically in terms of the massless GN model. 



continuum description mathematically equivalent to the 
symmetric (mo = 0) GN model. Dimerization, i.e., the 
fact that short and long bonds alternate, plays the role 
of (discrete) chiral symmetry breaking . Solitons and 
polarons are important for the physics of doping and 
electrical conductivity. They have originally been found 
by analogy with kink and kink-antikink baryons already 
known on the field theory side 12). By and large, the 
analysis of the phase diagram has evolved independently 
in the two fields and only recently merged into a consis- 
tent picture. 

If one starts to think about the massive GN model 
(mo 0), it is not too hard to identify its condensed mat- 
ter analogue: Conducting polymers with non-degenerate 
ground states like cis- PA (Fig. 2). Their theoretical de- 
scription was initiated by Brazovskii and Kirova [l3l | with 
the proposal that the gap parameter has two contribu- 
tions, a constant, "external" one arising from the basic 
structure of the polymer and an x-dependent, "internal" 
one due to electron-phonon coupling, 



A(x) = A e + Ai(x). 



(2) 
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If we identify A(x) with the scalar mean field S(x) and 
A e with the bare mass mo, we can immediately relate two 
problems from two different branches of physics. Baryons 
in the massive GN model correspond to polarons and 
bipolarons in the polymer case [f|. In the 80's and 90's, 
a lot of work was devoted to conducting polymers 0] 
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FIG. 2: Two inequivalent configurations of cis-PA, leading to 
a description in terms of the massive GN model. Configura- 
tion a) has lower energy than b). 

with the result that exact bipolaron lattice solutions were 
found by several authors. This strongly suggests that the 
massive GN model should exhibit a soliton crystal at fi- 
nite baryon density, much like the massless model. One 
has to be careful though, since in detail the use of the 
GN model is different in the two fields. In the relativistic 
quantum field case, one has to send the UV cutoff to infin- 
ity and the bare parameters Ng 2 and too to in a specific 
way, dictated by the renormalizability of the model. In 
polymer physics, the role of the UV cutoff is taken over 
by a band width W, a physical parameter (incidentally of 
the order of 10 eV). Likewise, the electron-phonon cou- 
pling A and the external gap parameter A e are physical 
observables to be taken from experiment. Nevertheless, 
in practice these differences do not seem to matter for 
many questions. We have noticed that condensed matter 
calculations are often done using approximations tanta- 
mount to the standard renormalization procedure in field 
theory. A similar remark applies to the role of N, the 
number of fermion species: Although N is as small as 2 
in polymer physics (the electron spin components), the 
adiabatic calculations done are in most cases indistin- 
guishable from the large N treatment on the field theory 
side. 

In the present work, we explore the possibility that the 
ground state of the massive GN model at finite baryon 
density has a crystal structure akin to the bipolaron crys- 
tal in polymers. We proceed as follows: We take advan- 
tage of the fact that the functional form of the periodic 
self-consistent scalar potential for non-degenerate poly- 
mers is known. It has been derived by various methods 
such as inverse spectral theory Il5l. Poisson summation 
of periodic sums of single polaro ns l l 6IIT7I fl8| or relation 
to the Toda lattice 0|. (Ref. [l!| contains a particu- 
larly thorough discussion and corrects misprints in some 
of the original papers.) We shall use the two-parameter 
functional form determined in these works as ansatz in a 
relativistic Hartree-Fock (HF) calculation, see Eqs. (|8I9|) 
below, but otherwise perform an independent calculation. 

Our main motivation is to generalize our recent calcu- 
lation of the phase diagram of the massless GN model 
to too 0. The most detailed investigation of the phase 
diagram of the massive GN model to date is the work of 
Barducci et al. 201. These authors assume that the chiral 



condensate is spatially constant and derive the phase di- 
agram for finite temperature and chemical potential. At 
T = in particular, they observe a 1st order phase tran- 
sition at a critical chemical potential depending on the 
bare fermion mass. Judging from what is known about 
the massless limit, these findings can probably not be 
trusted, since the existence of the baryon lattice has not 
yet been taken into account. 

Equipped with the structure of the single baryon in 
the massive GN model, the soliton crystal in the mass- 
less limit and a candidate for the self-consistent potential 
in the massive case borrowed from polymer physics, we 
should have all the ingredients necessary to clarify this 
problem. In the present work, we concentrate on finite 
baryon density at T = only. In Sec. II, we present the 
calculation of the ground state energy, assuming a peri- 
odic scalar potential of a definite functional form. Sec. Ill 
is devoted to the proof that the ansatz is indeed general 
enough to yield self-consistency, together with the deter- 
mination of its free parameters. In Sees. II and III, we 
rely on some results from Ref. Q where the same problem 
has been solved in the chiral limit, but nevertheless try 
to keep the presentation reasonably self-contained. Since 
the formulae are fairly complicated and their derivation 
is involved, we shall spend some time checking our results 
in various limits: The chiral limit too — > in Sec. IV, the 
low density limit in Sec. V and the high-density limit in 
Sec. VI. Once the formulae have passed all of these tests, 
we turn to exhibit their physics content and present some 
illustrative numerical results in Sec. VII. We summarize 
what we have learned in the concluding section, Sec. VIII. 
Two appendices arc devoted to the definitions of elliptic 
integrals and functions (Appendix A) and to the trans- 
lationally invariant solution of the massive GN model, 
needed here primarily to assess the stability of the soli- 
ton crystal (Appendix B). 



II. GROUND STATE ENERGY OF BARYONIC 
MATTER IN THE HARTREE-FOCK APPROACH 

Before describing in detail our calculations we give a 
short overview of what will be done next. In the present 
section, we shall compute the ground state energy density 
of baryonic matter in the relativistic HF approach. We 
assume right away a periodic ansatz for the scalar poten- 
tial S(x). This ansatz is borrowed from condensed matter 
physics and contains two free parameters, see Eqs. (JHEJ- 
It generalizes the scalar potential given in Eq. (6) of 
Ref. Q which was successful in the chiral limit of the 
GN model. The important feature of the ansatz which 
will enable us to carry out the calculation analytically is 
the fact that it yields the simplest, single gap Lame po- 
tential in the equivalent Schrodinger-type equation, just 
like in the chiral limit. Therefore many technicalities of 
the present calculation are unchanged as compared to 
this limit and a number of results from Ref. |3( can be 
taken over directly. The main difference between the two 
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calculations arises from the modified relation between the 
Dirac-HF single particle energies and the spectrum of the 
Lame equation. This renders summing up single parti- 
cle energies of occupied orbits somewhat tedious in the 
massive model. Although the final result of this section 
can still be given in closed analytical form, it requires 
incomplete elliptic integrals of all three kinds (see Ap- 
pendix A) as opposed to complete elliptic integrals in 
the chiral limit. This result will be presented in the form 
E g .s. = Ei +E 2 in Eqs. (|25I3U|) . where E\ is the sum over 
single particle energies of occupied states and E 2 the dou- 
ble counting correction to the interaction energy char- 
acteristic for the HF approach. It generalizes Eq. (24) 
in Ref. [3j to which it reduces in the chiral limit. The 
final expression for E g . s . of the present section still con- 
tains two undetermined parameters introduced through 
the ansatz for the scalar potential. The crucial question 
is then whether the ansatz is general enough to encom- 
pass the self-consistent potential. This will be answered 
positively in Sec. III. 

The first steps in the computation of the ground state 
energy density are the same as in the chiral limit. We set 
up the Dirac-HF equation which for the GN model Q 
reads 

1 b \^+ 1 S{x^ a {x)^u a ij a {x). (3) 

The scalar potential obeys the self-consistency equation 
which follows from minimizing the ground state energy 
SE g . s ./6S{x) = 0, 



1 



Ng- 



■(S(x) - rn ) 



occ. 

E 



1pa{x)lpa{x). 



(4) 



The sum runs over all occupied single particle states up 
to the Fermi energy including the filled negative energy 
states of the Dirac sea. In the present section, we concen- 
trate on evaluating the ground state energy for our ansatz 
for the scalar potential. The proof of the self-consistency 
condition which is technically the most difficult part 
of the whole calculation is deferred to Sec. III. 

From previous work |20j it is known that the HF prob- 
lem posed in Eqs. lol 111 has a solution with a uniform 
scalar potential S(x) = M (a physical fermion mass). 
For the convenience of the reader we have sketched this 
"standard" solution in Appendix B. There we also cover 
the vacuum case and show how to replace the dimen- 
sionless coupling g and the bare mass Too by physical 
parameters, thereby defining our renormalization proce- 
dure. The relation between bare and physical parameters 
[see Eqs. I|27I28I) ] will be crucial to cancel divergences 
from the Dirac sea, for instance on the right hand side 
of Eq. (0}. Apart from these well understood issues, the 
translationally invariant approach involves only straight- 
forward calculations, since all we have to do is to solve 
a free, massive Dirac equation. By contrast, a periodic 
ansatz for the scalar potential S{x) makes the HF prob- 
lem much more challenging. Fortunately we will see that 



owing to the specific form of S(x) the problem remains 
tractable. 

Returning to Eq. @, we use the following representa- 
tion of the 7-matrices, 

7° = -oi, 7 1 = iCT 3 , 7 5 = 7°7 1 = -a 2 . (5) 

In terms of upper and lower spinor components <j)± (sup- 
pressing the orbit label a from hereon), we obtain the 
pair of coupled equations 



which can be decoupled by squaring 

ox 1 ox J 
We now use the ansatz 

S(x) = AS(0 K = Ax) 



(0) 



(7) 



(8) 



where S is of the general functional form taken from work 
on bipolaron crystals in non-degenerate polymers [l5llT^ . 

EGam 

- _ sn£ + cnC+dn£+ + sn£_cnf_dn£_ 

o 27 27 ■ l y J 

sn z £\f — sn z £;_ 

Our conventions and shorthand notation for elliptic inte- 
grals and Jacobian elliptic functions are summarized in 
Appendix A, and 



e± = e ± s. 



(10) 



As is well known from quasi-one-dimensional condensed 
matter systems, the physics behind the appearance of pe- 
riodic structures is the Peierls instability — the system 
lowers its energy by generating dynamically an energy 
gap at the Fermi surface 21]. This was already observed 
in the chiral limit of the GN model where the state of low- 
est energy has a completely filled (for matter) or empty 
(for antimatter) valence band. Assuming that the gap is 
located at the Fermi surface, the scale factor A takes on 
the same value as in the chiral limit [K is the complete 
elliptic integral of the first kind, Eq. QA4J1 ]. 



A = 



2 Pf K 



(11) 



with the Fermi momentum pf related to the spatially 
averaged baryon density via p — Pf/ft- The two real pa- 
rameters which determine S(x) are the elliptic modulus 
k, suppressed in the notation of Eqs. (|9lllfl . and the shift 
5. With the ansatz 0, Eq. (JJJ) yields the simplest form 
of the Lame equation 



<9 2 

-— + 2 K 2 sn 2 e T 



b± = E(f>±, 



(12) 
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laying the ground for an analytical solution. The relation 
between Dirac- and Lame eigenvalues is 



where 77 is found to be 



±A^£ 



1 



sn 2 2<5 



(13) 



(14) 



The spectrum and eigenfunctions of the Lame equation 
(|12[1 are well-known [23. The energy bands k 2 < £ < 1 
and 1 + k 2 < £ are separated by a single gap, whereas 
the density of states inside the allowed bands is given by 

M 



dk 

d£ 



|E/K + k 2 -£ I 



2^(1 -£)(£- k 2 ){1 + k 2 - £) 



(15) 



Note that Eq. I|12|l is identical to what we had at mo = 0; 
at this point differences arise from the relation between 
Dirac single particle energies u> and Lame eigenvalues £, 
Eqs. ^13114(1 . and from reconstructing solutions of the 
original first order equations © where the shape of the 
potential S(x) enters. At a later stage, mo also appears 
in the self-consistency equation, see Eq. (QJ and Sec. III. 

In the present section, we focus on the ground state 
energy density 



E„ 



Ei + E 2 , 



(16) 



where E\ is the sum (or rather integral) over single par- 
ticle energies of all filled orbits (including the negative 
energy states in the Dirac sea) and E2 the HF double- 
counting correction. Proceeding like at mo = 0, we first 
sum over single particle energies, 



Ex = -2 A 2 



r £> d£ 


dk 


l £< 2^ 


d£ 



(17) 



with the integration limits (A/2 is the UV cutoff) 



l + K 

A 
2A 



2(1 -E/K) 



(18) 



[E is the complete elliptic integral of the second kind, 
Eq. (|A5|) ]. We consider antimatter here, filling only the 
negative energy states below the gap. This yields the 
same ground state energy as for matter, albeit in a sim- 
pler way. The integration in Eq. ((T7|) using (|T5)) is sig- 
nificantly more involved that in the massless case, but 
can be done along similar lines as in Ref . [l|j . The result 
is best expressed in terms of incomplete elliptic integrals 
F, E, LI of the first, second and third kind (see Appendix 



A) as follows, 
2tt 



A 2 



Ei = —a 



{£> 1) 



E 



- X [^~2 + K 2 )F( P ,q)+xE(p,q) 



k 2 / E 



— I 2— -2-»?)n(p,n,g), 



where we have introduced 



x = \A + v ■ 



dn 26 
sn2(5 



(19) 



(20) 



for a recurring expression to ease the notation. The ar- 
guments of the elliptic integrals are defined as follows, 



P = 




dn25 



(21) 



[k! = \/l — k 2 is the complementary elliptic modulus). 
Unlike in condensed matter physics, we are only in- 
terested in the limit A — > 00, so that we only need 
the asymptotic behavior of E\ for large £ > . The first 
term in Eq. (|19l) yields an irrelevant quadratic divergence 
—A 2 / (871") which can be dropped, as well as a finite term. 
In the 2nd and 3rd terms, we can replace p by p in E(p, q) 
and F{p,q), where 



p = —= = dn 2(5, 
\ n 



(22) 



all corrections being l/£> suppressed. The last term in 
Eq. (|19fl is more delicate, since II(p, n, q) has a logarith- 
mic singularity at p = p. Expanding 



1 



(23) 



v 7 ?! 2 v / n£ > 
and using a standard identity p4| in the form 

2(n - q 2 )(n - 1) 



n(p-e,n,g) 



1 



■ In: 



2 y (n — l)(n — q 2 ) ^Jn{n 2 — q 2 )e 
+ F(p,q)~u(pX,q] +0(e), (24) 



we finally arrive at 



2vr 



E 



El = 2--2 + ^--r? 



K 



2 g-2-,Vn, 



XE(P, q) 



+ j (2^-2-77) -n (/,(// r-</) 



\ [2^-2 + K 2 } Fip.a). 



(25) 
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The logarithmic divergence of the sum over single particle 
energies (the In A term in the 2nd line) will be cured once 
we add the double counting correction E 2 , to which we 
now turn. For finite bare quark mass mo it is given by 
(I = 2K is the spatial period in £), 



En 



I 



2Ng 2 I J 



d£(A§(0 - m ) 2 



(26) 



We invoke the vacuum gap equation [see Appendix B, 
Eq. 1B2II ] to eliminate the bare coupling constant, 



f 



f 



Ng 2 



(1 + m ) In A. 



(27) 



We use natural units in which the dynamical fermion 
mass in the vacuum is m = 1. In the massive GN model, 
there are two physical, renormalization group invariant 
parameters, the fermion mass m and 



7 = mo In A. 



(28) 



[Notice that our definition of 7 (called "confinement pa- 
rameter" in the condensed matter literature) differs by a 
factor of 7r from the definition in Refs. §,[2^. It turns out 
that the present definition simplifies the notation below.] 
We then find 



E2 = 2^ (7 + lnA) / d ^ 2(€) 



A7 



d£S(0- (29) 



Performing the integrations analytically, the final result 
for the double-counting correction becomes 



A 2 f F 



V 



A-f 



h-h 



II(K 2 sn 2 <5,K) 
K 



(30) 



with 



fx 
h 



cn 2 J + dn 2 <S + cn 2 <5dn 2 <5 

2 sn 5 cn S dn S 
2 cn 6 dn 6 



sn 5 



(31) 



[II is the complete elliptic integral of the 3rd kind, 
Eq. JUJ]. 

Eqs. (|25(l and (|30J) are the main result of this section. 
Upon adding up E\ and E 2 to get the ground state en- 
ergy density -E g . s ., the logarithmically divergent terms 
are cancelled and a finite result involving only physical 
parameters is obtained. It depends on the 4 parameters 
K, S,pf, 7. Out of these 4 parameters, p/ and 7 are deter- 
mined by the baryon density and the bare fermion mass, 
respectively, whereas n and 5 are so far unspecified pa- 
rameters of the trial potential. They will be determined 
in the following section by demanding self-consistency. 



III. SELF-CONSISTENCY CONDITION 

In the present section we check the self-consistency of 
the ansatz scalar potential Eqs. I|8I9|) . This ansatz con- 
tains two free parameters, K and 5. We proceed to show 
that the self-consistency condition (JIJ is satisfied pro- 
vided that k and 5 take on definite values depending on 
the density and the parameter 7. The main result of 
this section will be a pair of transcendental equations 
for k and S, Eqs. H46|) . which are equivalent to the self- 
consistency equation. 

The self-consistency condition for the scalar potential, 
Eq. Q), can be rewritten as 



-( 7 + lnA) + - = (32) 

IT IT ' ■ 



We have used Eqs. (|27|l and H28J1 to trade the bare param- 
eters Ng 2 and mo against 7 and the cutoff A. Our next 
task is therefore to evaluate the scalar density ijji/j for an 
arbitrary single particle solution of the Dirac-HF equa- 
tion and perform the sum over all occupied states. 
Since the first step involves only properties of the Lame 
eigenfunctions and is practically identical to the corre- 
sponding calculation in the appendix of Ref . Q , we skip 
some straightforward intermediate steps and just state 
the result, 



ijjij) 



An 2 



2w (dna - E/K) 



25 



(cn 2 a 



2 



The parameter a parametrizes the spectrum of the Lame 
Hamiltonian, see Eq. I|38() below. For the present pur- 
pose, it is actually more convenient to shift the variable 
^ by 5 and decompose the potential S, Eq. ©, into two 
terms, 



with 



51 = 

5 2 = 



S = S 1 + S 2 

s n $++cn£ ++ dn£ ++ 

sn 2 £ ++ — sn 2 £ 
sn £ cn £ dn £ 



(34) 



(35) 



sn 2 £++ - sn 2 £ ' 

We have introduced the doubly shifted argument 

£++=£ + 25. (36) 

The contribution to the condensate from a single orbit is 
then given by 



Using 



An 2 



•(dn 2 a-E/K) 



Ssn' 



dn 2 a = £ — k 2 



{I-E + k 2 



sn 2 a = H c 1 ^ 



(37) 



(38) 
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and summing over occupied states (as in Sec. II we con- 
sider antimatter) then yields 



r £ > d£ 


dfc 




d£ 



A_ 

"2% 



I 2 S + K 2 /i(£isn 2 £ + 5 2 sn 2 ^ 



(39) 



with the integrals 



h 



X 

2k 2 



h = \TL(p,n,q)-F(p,q)]. 

X 



(40) 



The arguments of the elliptic integrals are the same as 
above, see Eqs. Q21J1. Using once again the identity (|24() to 
isolate the logarithmic singularity of Tl(p, n, q) and taking 
the asymptotics with respect to A allows us to replace 



F(p,q) 
U(p, n, q) 



F(p,q), 

F(p,q)-U(p, ( K ') 2 ,q) 

1. X 1 A 
H — r In — — ■ . 

K Asj2 + Tj 



(41) 



Upon inserting the results I(39H41|) into the self- 
consistency condition (|32|l , the In A term drops out once 
again and we arrive at the finite equation 



XI 
An 2 



(5-1) = F{p, q) Sisn^ + S 2 m^+ 



U(p,(K') 2 ,q) + ^ln(A^2T^) 



(42) 



S. 



To clarify the content of this condition, we re-arrange it 
slightly and introduce coefficients C\ : 2 as follows, 



S- 



c 



7 = Ci 

d = Ahi{A^2 + n)+A-/ 



-± (S lS n 2 £ + S 2 sn 2 £++ 
An 2 



n(p, (K r ) 2 ,q) 



C 2 = F(p,q). 

X 



(43) 



At first glance, it seems unlikely that Eq. (|43|l can be 
solved, since a constant term and two different functions 
of x appear. However, these functions are not linearly 
independent. Using addition theorems for cnoidal func- 
tions, one can verify the following identity, 



S - K 2 sn 2 2<5 (£isn 2 £ + S 2 sn 2 ^ 



cn2<5dn2<5 
sn2(5 



(44) 



Comparing Eqs. (|43|l and 144(1 . the self-consistency condi- 
tion may be turned into the following two x-independent 
conditions, 



Ci 
C2 
Ci 



cn2 ( 5dn2 ( 5 
sn2<5 ' 

-K 2 sn 2 25, 



(45) 



with Cx,2 as defined in Eq. (|43|l . For 7 7^ 0, these equa- 
tions can be cast into the somewhat simpler form 



= Acn2SF(p,q) - 7sn 2 26, 



(46) 



7 



Acn2S 



K 2 n (p, (k') 2 , q)+ X (l + H^V^+V)) 



They determine k and S for given pf and 7. 

After carefully identifying all the variables, the self- 
consistency equations (|4^)l agree exactly with Eqs. (5,6) 
in |l9j |. confirming the 1:1 mapping from the theory of 
non-degenerate conducting polymers to the massive GN 
model. If we can solve the pair of transcendental equa- 
tions (|46jl . we have found a solution of the HF equa- 
tions and hence a candidate for the ground state of bary- 
onic matter. Another solution of the HF equations is 
the translationally invariant one discussed in detail in 
Ref. [23] and summarized in Appendix B. Which solu- 
tion is favored is then simply a question of the energy 
which can be computed from Eqs. (|25I3U[) of Sec. II. The 
numerical results will be presented in Sec. VII. The fol- 
lowing three sections are devoted to testing the general 
formalism in simple special cases where we know the an- 
swer from other sources. 



IV. CHIRAL LIMIT 

As a first and most trivial test, let us check the above 
formalism in the chiral limit (7 = 0) against our previous 
work. Consider the self-consistency equation, Eq. I|45[) . 
The first equation can be solved by 



6= 2> 



cn28 = 0. 



Consequently, 

P = X = k', q=l, V 



(47) 



(48) 



and the elliptic integrals are reduced to elementary func- 
tions, 



F(k', 1) = artanhft', 
H(k',( K ') 2 ,1) 



(49) 



If,, K 2 , 1 + K 1 

— 77 k in ; — -77 + in 

2k 2 V 1 + («') 



1 - k' 



The second equation of Eqs. (|45|l then yields simply 

An = 1, (50) 

relating k and pf. Likewise, the ground state energy 
Eqs. H25I30I) simplifies tremendously in the limit 7^0, 



4tt I K 



E = — I 4— - 2 + K 2 



A 2 ( E \ 
— ( 2— - 2 + k 2 I ln^K. 
2n \ K J 

(51) 



Finally, starting from Eqs. (|34I35|1 and using standard 
identities for cnoidal functions of the shifted argument 
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£ + K, the self-consistent potential can be manipulated 
into the form 



5' 



2 sn£cn£ 
v dn£ ■ 



(52) 



All of these formulae agree with the known results for 
7 = 3]. A comparison of the full equations of Sees. II 
and III with the reduced ones in the present section 
shows a significant (but apparently unavoidable) increase 
in complexity due to the finite bare fermion mass. 



V. LOW DENSITY LIMIT AND SINGLE 
BARYON 

A second obvious way of testing the formalism is the 
low density limit, where we expect to recover the prop- 
erties of single baryon s in the massive GN model. We 
have to expand re = y/l — (k 1 ) 2 around k 1 — 0. To lead- 
ing order, we may neglect all power corrections in k' but 
must keep the logarithmic singularity ~ In k' whenever 
it is present. In this way, we find 



p w sech2<5, 



0, 



(53) 



and hence all incomplete elliptic integrals reduce to the 
same elementary function, 



{F(p,q),E(p,q) 1 U(p, (n') 2 ,q)} -> arcsinp. 
The self-consistency conditions (|46|) become 



= 2pfp&rcsmp\n 7T7tanh 25 



(54) 



(55) 



2p/P hi- 
re 



with the definition 

T = In 



1 



arcsmp - 



sinh 28 



( 7 + r) 



2p/ln(4/«'; 
7r tanh 26 



(56) 



These equations can be solved parametrically for 5 and 
re' as follows, 



1 



artanhy, 



2 

4 exp 



2p//' 



(57) 



where 



sin6>, 

7T 



9 — 7 tan 0, 



(58) 



in perfect agreement with the corresponding equations 
for the single baryon (with completely filled valence level) 
in the massive GN model 



Using our full expression for the energies E\ and E 2l 
Eqs. (|25I30|I . the limit re' — > yields the finite parts 



£1 



\_ pi ( 2yi 2y 

ATT IT \ 7T 7T 



E 2 = -f + -/f^ + -I ar tanh y 
2ir n \ n n 



(59) 



The logarithmically divergent pieces have been omitted, 
since they anyway cancel in the sum. The total energy 
density can then be represented as 



E, 



g.s. 



Ei+E 2 



E v 



P M E 



(60) 



where the vacuum energy density agrees with Eq. (|B4|) 
in the appendix, p = pf/w is the mean baryon density 
and Mb is the known baryon mass 5] . This is indeed the 
expected low density behavior for widely spaced baryons 
and an important additional test of the formalism. 

Finally, we check the shape of the scalar potential in 
the low density limit. Taking re — » 1 in Eq. Q, we find 



lim S(0 = coth2(T 



sinh 28 



t->i 



cosh 2 £ + cosh 2 5 — 1 



(61) 



Let us compare this expression to the profile of the single 
baryon at m ^ 0, 



Sb(x) = 1 + y [tanh(ya; — cq) — tanh(ya; + cq)] 
1 



— artanht/. 
2 y 



(62) 



yx, 



Co, 



(63) 



Co 

If we identify 

A = y, £ 

this can be replaced by 

S B {Q = coth 25 + tanh(£ -5)- tanh(£ + 5) (64) 

in agreement with Eq. H61|l. One can easily convince 
oneself that the relations I|62I63|) are equivalent to the 
self-consistency equations (1571) . 



VI. HIGH DENSITY LIMIT AND 
PERTURBATION THEORY 

In the high density limit pf — > 00, one would expect 
that all interaction effects can be treated perturbatively. 
As is well known, band formation in a periodic potential 
requires almost degenerate perturbation theory (ADPT). 
Here we closely follow a similar calculation carried out by 
us in the context of the massless GN model 0- There 
are two differences: In the double-counting correction, 
we have to take into account the bare quark mass, and 
we allow for £0 7^ in addition to S±± 7^ [Se are the 
Fourier components of the periodic potential S(x)]. Thus 
our present ansatz for S(x) is 



S(x) = S + 5ic i2p ^ + SUe 



-i2p j-x 



(65) 
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The potential has a spatial period a equal to the inter- 
baryon distance, i.e., the inverse density 



7T 

Pf' 



(66) 



This is the reason why the lowest non-vanishing mo- 
mentum which appears in the Fourier expansion (|65|l is 
lix I a = 2pf. Without loss of generality, we may assume 
that So and S\ = S-% are real (a phase in Si corresponds 
to a translation). The sum over single particle energies 
(for the antimatter case) in 2nd order perturbation the- 
ory (PT) is given by 



Ei 



A/2 



dk 



2k 



SI 



2(k+p f ) 
yJ(k- Pf )* + S* + (k- Pf )} 

si 



A 2 
1 

Til 



Pf 



2lT 47T 



1 

2^ 



*o 2 



2Sf) In ■ 



A 



Si 



S^ln^ + ^lnSi. 



(67) 



We had to invoke ADPT only for the term which would 
blow up in naive 2nd order PT, 



si 



, {L ._ i)j] y (k-p f y + Sf + (k-p f ) (68) 

This "recipe" has been derived in Eq. (2.11) of Ref. Q 
by simply comparing 2nd order degenerate and non- 
degenerate perturbation theory. In Ref. [2J, the states 
which are almost degenerate in a periodic potential 
are explained in more detail, following the standard 
weak binding approximation from solid state physics (see 
e.g. [1^|). The double-counting correction for the poten- 
tial (|d3|| becomes 



E 2 



1 1 



dx(S(x) 



2Ng 2 L 

^{Sl + 2Sf)(7 + In A) 



mo) 2 

7 



So- 



(69) 



Adding Eqs. (|67|l and 169J1 then yields the approximate 
energy density 



A2 
'8tt 

si 

2ir 



2tt 
;<2 



_L 

2tt 



S — 25q) 



^(ApfSJ + lsl 



Minimizing with respect to 5*0 and Si, we find 

SoK 2 P/)+7]-7 = 0, 
5i [2 1 + \n(4 Pf S 1 )}=0. 

The first equation has the unique solution 

o _ 7 



7 + ln(2p/) 



(70) 



(71) 



(72) 



-0.002 



-0.004 



AE 



-0.006 



-0.008 



-0.01 



-0.012 




FIG. 3: Difference in energy density between crystal and 
translationally invariant Fermi gas. From top to bottom: 
7 = 2.3,0.75,0.3,0.1,0.01. Dashed curves: asymptotic be- 
havior according to Eqs. 1601 and I74H . for 7 = 0.01. 



in agreement with the leading term in Eq. Ijf312jl of the 
appendix. The 2nd equation has two solutions: S\ = 
0, corresponding to unbroken translational invariance as 
discussed in Appendix B, and 



Ap f 



-27 



(73) 



for the soliton crystal. Comparing the energy densities 
of these two solutions, 



E s . s .(Si^0)-E g . s .{Si=0) 



1 



64tt^ 



-47 



(74) 



we learn that the crystal is favored, but the energy differ- 
ence decreases rapidly with increasing 7. Eqs. (|73|l and 
(|74J) agree with our previous results if, in addition to the 
high density limit, we take the chiral limit 7 — > 0. 



VII. NUMERICAL RESULTS 

Before computing any observable, we have to solve the 
self-consistency equations (|45|) . Choosing the mean den- 
sity and the bare fermion mass, or, equivalently, p/ and 
7, these two transcendental equations yield the two un- 
known parameters k (elliptic modulus) and S (shift pa- 
rameter) of the trial potential defined in Eqs. (jST lllJI . 
Eqs. (|45|l always have two different solutions, a trans- 
lationally invariant one (k = 1) and a spatially modu- 
lated one (0 < k < 1). In the case of the translationally 
invariant solution, there is an additional complication al- 
ready familiar from the massless limit: At low densities, a 
mixed phase appears, characteristic of a first order phase 
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transition. Therefore, one cannot simply take the HF 
energy at face value, but has to introduce a second vari- 
ational parameter describing the amount of space filled 
with "droplets" of baryonic matter. Alternatively, one 
can work with a chemical potential and use the grand 
canonical potential at T = as in Ref. [20|. leading to 
the same results. More details are given in Appendix B. 
The energy density of the crystal phase can be computed 
by inserting the self-consistent values for k and d into 
Eqs. (|25I3(J|) . We then find that the crystal is energeti- 
cally favored at all densities and all bare quark masses. 
In the low and high density limits, this can be shown 
analytically, see Eqs. (|B23IB24|1 for the limit pj — ► 
and Eq. Ij74(l for pf — > oo. In between, we have to 
compute the energy difference numerically (where "nu- 
merically" in this context simply means using floating 
point commands in Maple, thus getting any desired ac- 
curacy). Some illustrative results are shown in Fig. 3 for 
five different values of 7. Together with the lowest curve 
(corresponding to 7 = 0.01), we have also plotted two 
dashed curves corresponding to the (analytical) asymp- 
totic behavior at small and large densities. We observe 
that the agreement of the full calculation with the asymp- 
totic curves is excellent, with a narrow crossover region 
where the curve changes rather abruptly from ~ pf to 
~ pj 2 behavior. 

On the basis of such calculations, we conclude that 
the lattice solution is stable at all densities and quark 
masses. In polymer physics, the situation is somewhat 
different, and stability of the bipolaron lattice has been 
found in certain regions of parameter space only |19| . 
The difference can be traced back to the fact that in 
the condensed matter case, the cutoff A is not sent to 
infinity, but is a physical parameter (the bandwidth of 
the undimerized polymer) . Here the GN model is used as 
an "effective theory" , like the Nambu-Jona-Lasinio model 
in particle physics. One therefore has more parameters 
to vary, see e.g. Fig. 4 of Ref. 0|- The horizontal axis 
in this plot is labelled A e /A, corresponding to rrio/m 
in our case. Since the bare fermion mass goes to in 
the limit A — > 00, we can only compare our calculation 
with the point A e /A = 0, where there is indeed stability 
also in the polymer case. The unstable region for finite 
values of mo is simply not accessible if one treats the GN 
model strictly as a renormalizable relativistic quantum 
field theory. 

Having established the stability of the crystal solution, 
let us now illustrate how a finite bare quark mass modi- 
fies the self-consistent scalar potential. This is exhibited 
in Figs. 4 and 5 at low and high density, respectively. 
The deepest curves always correspond to 7 = 0, where 
the potential oscillates symmetrically around zero. This 
is actually a remnant of the original discrete chiral sym- 
metry of the model. Translational invariance and the 7 
transformation both break down, leaving an unbroken 
discrete symmetry (translation by half a period, com- 
bined with a 7 5 transformation). Since the massive GN 
model does not have discrete chiral symmetry in the first 




x 



FIG. 4: Self-consistent scalar potential S(x) versus x for p/ = 
0.3. From bottom to top: 7 = 0, 0.01, 0.1, 0.3, 0.75, 2.3. 



0.6 ~- 



0.4 - 




-2 -1.5 -1 -0.5 0.5 1 1.5 2 



FIG. 5: Same as Fig. 4, but for pf = 2.5. Note the different 
scale on the z-axis. 



place, one would not expect the same behavior here. In- 
deed, the potentials now wiggle around a finite value, 
with a less symmetric shape. In the massive case, contin- 
uous translational invariance is broken down to a resid- 
ual discrete translational invariance. As we turn on the 
symmetry violation parameter 7, the potential oscillates 
with decreasing amplitude around a value close to the 
mass M which the fcrmions would acquire in the transla- 
tionally invariant solution, eventually leaving only a very 
weak modulation of a large scalar potential in the heavy 
fermion limit. It is surprising that such a variety of po- 
tential shapes in the Dirac equation can all be reduced 
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FIG. 6: Baryon density versus chemical potential for trans- 
lationally invariant solution, showing a first order transition. 
From left to right: 7 = 0.01, 0.1, 0.3, 0.75, 2.3. 



FIG. 7: Same as Fig. 6, but for crystal solution. Here, a 
continuous, second order phase transition at fi c = Mb is seen 
for all values of 7. 



to the standard single gap Lame equation. 

The last result which we should like to show is how the 
density varies with the chemical potential. The chemical 
potential at T = can be obtained by differentiating the 
energy density with respect to the mean fermion density, 



8E. 



A* = 



dp 



(75) 



If we assume unbroken translational invariance (Fig. 6), 
we find discontinuities in these curves, confirming the 
result of Ref. [2(j about a first order phase transition. 
Repeating the same calculation for the crystal solution 
(which is the stable one), all the curves become continu- 
ous, signalling a 2nd order phase transition (Fig. 7). The 
critical chemical potential in this latter case coincides 
with the baryon mass, as expected on general grounds. 
By contrast, the first order transition in Fig. 6 happens 
at a chemical potential which has at best the meaning of 
an approximate baryon mass in a kind of droplet model, 
cf. Appendix B. 



VIII. SUMMARY AND CONCLUSIONS 

In this paper, we have derived analytically the ground 
state of cold, dense matter in the massive GN model. 
This task was greatly simplified by the lucky circum- 
stance that the problem is mathematically closely re- 
lated to another problem which had already been solved 
in condensed matter physics, namely the bipolaron lat- 
tice in non-degenerate conducting polymers. Since we 
could not be sure that subtle differences between the 



polymer problem and relativistic quantum field theory 
don't matter, we only took over the functional form 
of the self-consistent potential from literature and car- 
ried out an independent HF calculation with this two- 
parameter ansatz. Our results fully confirm the 1:1 map- 
ping between two seemingly unrelated physics problems. 
Whereas a similar relationship had been noticed before 
by several authors in the context of degenerate polymers 
and the massless GN model, this time we were able to 
take full advantage of the correspondence in a situation 
where we could never have guessed the correct potential. 

As compared to our previous work on the massless GN 
model, the calculations become significantly more com- 
plex once one turns on the fermion mass parameter 7. In 
order to write down the ground state energy, one now 
needs all 6 standard complete and incomplete elliptic 
integrals. The results exhibit a wider variety of scalar 
potentials oscillating not around 0, but around a finite 
value. Since the Dirac equation with all of these poten- 
tials can still be reduced to the single gap Lame equation, 
the calculations remain tractable. The crystal is energet- 
ically favored over the translationally invariant solution 
at all densities and quark masses as a result of the Peierls 
effect with gap formation at the Fermi surface. Like in 
the chiral limit, the alleged first order phase transition 
discussed in the literature disappears in favor of a 2nd 
order phase transition. The critical chemical potential 
now coincides with the baryon mass, an important theo- 
retical constraint which had been violated by the trans- 
lationally invariant solution. It would be very interesting 
to extend this work to finite temperature and determine 
the full phase diagram of the massive GN model. What 
role does the soliton crystal play at finite temperature, 
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and what happens to the first order transition line ending 
at a critical point |20j? Provided that the ansatz for the 
scalar potential is flexible enough to describe the finite 
temperature case as welL this would be rather straight- 
forward, following Ref. |4| for the massless limit. As an 
additional bonus, such a calculation might even be of 
some interest for polymer physics. 

Having seen the quantitative correspondence between 
polymers and the GN model, one wonders whether other 
results could be traded between the Phys. Rev. D and 
Phys. Rev. B communities and might shed some new 
light on the respective physics questions. One such issue 
might be asymptotic freedom, the property to which the 
GN model owes its existence in the first place. Is there 
anything analogous to signals of asymptotic freedom (like 
scaling in deep inelastic scattering) in polymer physics, 
perhaps in optical properties of conducting polymers? 
Conversely, in the condensed matter literature, the sym- 
metry breaking parameter 7 is known as "confinement 
parameter". Since the ground state is no longer degen- 
erate, the effective kink-antikink potential rises linearly 
with distance. Could this mean that the GN model can 
teach us anything about confinement in quantum chro- 
modynamics? It will be interesting to follow further such 
speculations and identify questions where an intensified 
dialogue between particle and condensed matter physi- 
cists might be fruitful. 

We should like to thank Oliver Schnetz for mathemat- 
ical support and useful conversations. We also thank the 
referee of Physical Review D for helpful suggestions of 
how to improve the readability of an earlier version of 
this paper. 



APPENDIX A: ELLIPTIC INTEGRALS AND 
JACOBI ELLIPTIC FUNCTIONS 

Since different ways of writing down elliptic integrals 
and Jacobi elliptic functions arc in use, we briefly sum- 
marize our conventions. As a rule, we use the Legendre 
normal form for elliptic integrals. Our notation of the 
arguments is the same as in Maple. 

• Incomplete elliptic integral of the first kind 

1 



F(z,k) 



dt- 



lo VT^Vl - n 2 t 2 
• Incomplete elliptic integral of the second kind 



(Al) 



• Incomplete elliptic integral of the third kind 

1 



II(z, n, k) = / dt — 
1 ' I (1 



(1 - rrf 2 )vT^Vi - K 2 t 2 



(A2) 



(A3) 



• Complete elliptic integral of the first kind 

K = K(«) = F(l, k) (A4) 

• Complete elliptic integral of the second kind 

E = E(k) = E(l, k) (A5) 

• Complete elliptic integral of the third kind 

n(n, k) = 11(1, n, k) (A6) 

In Jacobi elliptic functions, we suppress the elliptic mod- 
ulus k throughout this paper, 

snu = sn(u,ft), cnu = cn(u, k), dnu = dn(w, k). (A7) 

Following the convention implemented in Maple, we have 
denoted the second argument by k rather than k 2 . 



APPENDIX B: TRANSLATIONALLY 
INVARIANT APPROACH 

For the sake of comparison with the crystal, we also 
need results for the translational invariant HF solution 
pcj . Let us start with the vacuum and recall the basic 
equations 0. The vacuum energy density reads 



A/2 



E = -2 / ^Vp 2 +m 2 + 



2tt 



(m — m^Y 
2Ng 2 



A 2 m 2 m m (to — mo) 2 
= ~8^~ 2^ ln A + 2Ng 2 ' (B1) 

Varying with respect to the Fermion mass to yields the 
gap equation 



7T , A 

7 + hi — 



with 



Ng 



_ toq A 

TO TO 



(B2) 



(B3) 



From now on, we set to = 1 as in the main text. For the 
vacuum energy density at the minimum, we get 



E = 



A 2 1 2_ 

87T 47T 27T ' 



(B4) 



The energy density of the massive Fermi gas differs from 
Eq. I|B1J| only in the integration limits, 



e = -2 r ^vp^ + 2Ng2 



(M - too) 2 



ps 



A 2 M 2 M 2 M (M-to ) 2 
" 8^ ~ 4tT + ~2T 11 T + 2Ng 2 



2?r 2tt 



(B5) 
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with 



V) 



M 2 



(B6) 



Eliminating the bare coupling constant with the help of 
the vacuum gap equation, the logarithmic divergence dis- 
appears and we are left with 



8tt 



M 2 

47T 



+ 



M 2 

2vr 



M 

— J 7" 



-^+_ln(p / + e/ ). 

(B7) 

M can be obtained by minimizing E with respect to M, 

(M - 1)7 + M In ( P/ + e/ ) = 0. (B8) 

The difference between energy density and vacuum en- 
ergy density is finite, 



E-E 



(M 2 - 1) (M - 1) 



47T 



2tt 



-7' 



2tt 



(B9) 



This last equation is only valid at the minimum, since we 
have made use of Eq. (|B8|) to simplify the expression. 

In general, Eq. IjBSjl can only be solved numerically. In 
some limiting cases, one can solve it by a series expansion 
for M and compute the corresponding energy. Consider 
the following limits: 



• 7 - 



at fixed pf (heavy quark limit), setting 



M 



1 hx{p f + e s ) 

7 



r 



(BIO) 



Energy density (without vacuum subtraction), 



E ps — ! - 



Pf e f 



1 

47T 2 ,T 

. 2, 



— ln(p f + e f ) 



• Pf ~ 

M 



00 at fixed 7 (high density limit): 



1 



7 + ln2p/ Ap) (7 + 1H2P/) 4 ' 



(Bll) 



(B12) 



P/ 1 7 2 1 7 4 ffil3l 

2tt 2tt 7 + In 2p f l&irpj, (7 + In 2p/) 4 ' 1 ^ 



As described by Barducci et al. |2(j, the translationally 
invariant solution undergoes a first order phase transi- 
tion at a certain critical chemical potential or density. 
In the HF solution described above, this manifests itself 
through the fact that the energy- versus-p/ curve starts 
out as concave at pj = and becomes convex only above 



a certain critical Fermi momentum p c j. Before drawing 
conclusions about the stability of the crystal, we have to 
take this fact into account. 

Just like at 7 = 0, we make a mixed phase variational 
ansatz 9]: A fraction A of space contains all the extra 
fermions ( "droplets" ) , which therefore have an increased 
Fermi momentum 



r/ - EL 

f ~ A ' 



(B14) 



The energy density, subtracting the vacuum contribution, 
is 



AE 



X 



r m 2 


(M 2 






V 2tt 


-t) 



(B15) 



P'f e 'f M 2 , / / M 1 7 1 
2tt 2tt v 1 47r 2n J 



with 



e' f = .^M 2 + (p' f y 



Vary with respect to A, 
dAE 



dX 



= ^ AE - 



■p / ^- = 0. 



dp* 



(Bi6) 



(B17) 



/ 



The right hand side can be interpreted as construction of 
the convex hull of the curve energy density versus fermion 
density (remember the geometrical meaning of the Leg- 
endre transform). The solution p'^ of this equation at 
A = 1 is the critical Fermi momentum. More explicitly, 
Eq. llBT7|l reads 



(B18) 
(B19) 



Eqs. I|B18(I and i|B19(l determine the mixed phase. To 
solve them, proceed as follows: Eliminate the ln-term 
in Eq. (|B18|) with the help of Eq. (|B19|) and solve the 
resulting equation for p'^ , 



(l+Af+2 7 -27M)(l-Af)-2p' / e' / +2M 2 In (p' f + e' f ) = 
Vary with respect to M [cf. Eq. JBS}], 

(M-l)7 + Mln(p' J + 4) =0. 



-M 2 + ^M 4 + (1 - M) 2 (l + M + 2 7 ) 2 

2 ' 

(B20) 



(P'f) 2 



Insert this result into Eq. (|B19|) and arrive at an equa- 
tion relating M and 7. This equation always has a non- 
trivial solution, which has to be determined numerically. 
Analytical results can be obtained in limiting cases of 
interest, 



M 



M 



:7 



In 2 ' (In 2 
for 7 — > 0, 
3 3 



2 4[101n2 + (ln2) 2 + 2] 3 
^7 H ru~^T2 7 



1 



G3 



2 7 2 7 3 4Q 7 4 



(In 2) 4 
for 7 — > 00. 



(B21) 
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By inserting these expressions into Eqs. (|B20|) . (|B15(1 . we 

can determine the behavior of the energy in the mixed 
phase, as well as the critical point. Find for small 7 



Pf f 1 



1 



AE w £^ _ + _ 7 - 



4 + In 2 



7T 1 V2 V2 2\/21n2 



■41n2+(ln2) J 
2\/2(ln2) 2 



7 2 , (B22) 



and for large 7 



AE - -I^A+A 

7r [ 87 2 47 3 
3 3 69 
2^ ~ 27^ ~ 80 7 3 ' 



(B23) 



Pf 



The expressions in curly brackets in Eqs. I|B22|) and (|B23 



are the critical chemical potentials (or, equivalently, the 
baryon masses in a kind of "bag model" |9J). Compare 
this with the true baryon masses in the two limits [cf. 
Eqs. , 



M B 
M b 



1 i I 4/ 



1 - 



24 7 2 !2 7 3 



for 7 
for 7 — 



(B24) 



In these limits, the baryon mass is below the bag model 
mass. One can easily check numerically that this is in fact 
true for arbitrary values of 7. This proves that the trans- 
lationally invariant solution is always unstable against 
formation of a kink-antikink crystal. 
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